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Summary 

Four turbulence models are described and evaluated for transonic flows over the High-Speed Research/ 
industry baseline configuration known as Reference H by using the thin-layer, upwind, Navier-Stokes solver 
known as CFL3D. The turbulence models studied are the equilibrium model of Baldwin-Lomax (B-L) with 
the Degani-Schiff (D-S) modifications, the one-equation Baldwin-Barth (B-B) model, the one-equation 
Spalart-Allmaras (S-A) model, and Menter’s two-equation Shear Stress Transport (SST) model. The flow 
conditions, which correspond to tests performed in the National Transonic Facility (NTF) at Langley 
Research Center, are a Mach number of 0.90 and a Reynolds number of 30 x 10 6 based on mean aerody- 
namic chord for angles of attack of 1°, 5°, and 10°. The effects of grid topology and the representation of the 
actual wind tunnel model geometry are also investigated. Computed forces and surface pressures compare 
reasonably well with the experimental data for all four turbulence models. 

Introduction 

Increasing computer capacity and the development of efficient numerical methods allow the solution of 
the Navier-Stokes equations for practical aerodynamic applications such as a next generation supersonic 
transport. Of crucial importance for the accuracy of the results is the performance of the turbulence model 
because the viability of a supersonic transport is, in part, based upon the ability to predict drag to within one 
count. It has been observed during the High-Speed Research (HSR) program that computational solutions 
have been inconsistent, showing dependence on turbulence model and model orientation, among other fac- 
tors (ref. 1). With the increasing complexity of the flows under investigation, the demand for reliable turbu- 
lence models is as high as ever; these models must perform well for attached flow and at least the onset of 
separated flow that occurs over an angle-of-attack range. The search for a new or improved turbulence model 
that accurately predicts both attached and separated three-dimensional flow fields is complicated by inherent 
limitations in the computational fluid dynamics (CFD) codes such as truncation error, grid density, the type 
of differencing scheme employed, and for central-difference schemes, the amount of artificial dissipation 
added for numerical stability. 

This paper investigates the ability of four different turbulence models to predict accurately the transonic 
flow over the High-Speed Research/industry baseline configuration known as Reference H. This configura- 
tion is an attached flow airplane that goes to separated flow at high angles of attack. The thin-layer Navier- 
Stokes solver CFL3D (ref. 2) was run using the Spalart-Allmaras (S-A) (ref. 3), the Baldwin-Lomax (B-L) 
(ref. 4) with the Degani-Schiff (D-S) modifications (ref. 5), the Baldwin-Barth (B-B) (ref. 6), and Menter’s 
Shear Stress Transport (SST) turbulence models (ref. 7) for this study. Each model is discussed in detail 
herein. No attempt was made to tune these turbulence models for this configuration. 


Symbols 


a l(SST) 

constant used in Menter’s Shear Stress Transport turbulence model 

b 

model span, in. 

C D 

drag coefficient 

C Dp 

pressure component of drag coefficient 

Q>v 

viscous component of drag coefficient 

C L 

lift coefficient 

Cm 

pitching-moment coefficient 

C P 

pressure coefficient 
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mean aerodynamic chord, in. 
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local chord, in. 

function in Baldwin-Lomax turbulence model 

distance closest to any wall in Spalart-Allmaras turbulence model, in. 

total energy per unit volume, ft-lb/ft 3 

function in Baldwin-Lomax turbulence model 

inviscid fluxes 

functions used in Menter’s Shear Stress Transport turbulence model 

viscous sublayer damping terms used in Baldwin-Barth turbulence model 

functions in Spalart-Allmaras turbulence model 

viscous fluxes 
Jacobian matrix 
kinetic energy, slug-ft 2 /sec 2 
lift-to-drag ratio 
free-stream Mach number 

production term used in Baldwin-Barth turbulence model 
conservation variables 
dynamic pressure, psf 

Reynolds number based on mean aerodynamic chord 
turbulence Reynolds number 

working variable related to Rer, Re r = (Re r ) A (Re r ) in 
Baldwin-Barth turbulence model 
time, sec 

Cartesian velocities in x, y, and z directions, respectively, ft/sec 

velocity vector used in Menter’s Shear Stress Transport turbulence model, 
ft/sec 

Cartesian coordinates, in. 

position vector used in Menter’s Shear Stress Transport turbulence model, 
in. 


distance to wall or wake cut in Baldwin-Lomax turbulence model 


y* 

a 

P(SST)’ P(SST) 

r 

Y(SST) 


law-of-the-wall coordinate, 





angle of attack, deg 

constants used in Menter’s Shear Stress Transport turbulence model 
Klebanoff intermittency factor used in Baldwin-Lomax turbulence model 
constant used in Menter’s Shear Stress Transport turbulence model 
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p 

molecular viscosity, slugs/ft-sec 

Hr 

molecular eddy viscosity, slugs/ft-sec 

Hw 

molecular viscosity at wall in Baldwin-Lomax turbulence model, 
slugs/ft-sec 

V 

kinematic viscosity, ft 2 /sec 

V (S-A) 

working variable in Spalart-Allmaras turbulence model 

v r 

kinematic eddy viscosity, ft 2 /sec 

5.T1.C 

general curvilinear coordinates 

P 

density, slugs/ft 3 

P*v 

density at the wall in Baldwin-Lomax turbulence model, slugs/ft 3 

^eCB-B) 

function used in Bald win- Barth turbulence model 

°*<SST)’ °co(SST) 0 cu2(SST) 

constants used in Menter’s Shear Stress Transport turbulence model 

x y(SST) 

Reynolds stress tensor used in Menter’s Shear Stress Transport 
turbulence model 


local shear stress in Baldwin-Lomax turbulence model, lb/in 2 

n 

magnitude of vonicity. J ( g - + (| - |£f + (£ - 1) 2 

© 
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Abbreviations: 

specific dissipation rate, 1/sec 

AF 

approximate factorization 

B-B 

Baldwin-Barth 

B-L 

Baldwin-Lomax with Degani-Schiff modifications 

C-O 

C in streamwise direction; 0 in spanwise direction 

CPU 

central processing unit 

D-S 

Degani-Schiff 

Exp. 

experiment 

FAS 

full-approximation scheme 

HSR 

high-speed research 

LE 

leading edge 

MW 

megawords 

NAS 

Numerical Aerospace Simulator 

o-c 

0 in streamwise direction; C in spanwise direction 

PDE 

partial differential equation 

resid 

residual 

S-A 

Spalart-Allmaras 

SST 

Shear Stress Transport 



Subscripts: 

i 

y 

max 

min 

o 

t 

x,y 

w 

wake 


inner region values 
tensor notation indices 
maximum value 
minimum value 
outer region values 
differentiation in time 

differentiation in x and y directions, respectively 
conditions at wall 
conditions in wake 


Experimental Study 

Tests of the High-Speed Research/industry baseline configuration known as Reference H have been con- 
ducted in the National Transonic Facility (NIT) at Langley Research Center. The 2.2-percent scale model of 
the clean wing-body configuration is shown in the NTF test section in figure 1 ; the model has a wing span of 
34.23 in., a mean aerodynamic chord (c ) of 22.71 in., a reference area of 3.436 ft 2 , and pitching-moment 
data have been referenced to the 50-percent c location (x = 46.445 in.). The leading-edge sweep of the wing 
is 76° inboard of the 23-percent semispan, 68.5° from the 23-percent to the 52-percent semispan, and 48° 
outboard of the 52-percent semispan. The leading edge is blunt inboard of the 52-percent semispan and sharp 
outboard of the 52-percent semispan. Figure 2 shows the location of the six spanwise pressure rows on the 
wing, and figure 3 shows the location of the three chordwise pressure rows on the wing of this geometry. 
Although the shaded area in these two figures was part of the metric model in the wind tunnel experiment, it 
was not part of the computational model in order to make the grid generation easier. 

The NTF is a cryogenic, pressure wind tunnel with independent control of pressure, temperature, and fan 
speed. The independent controls allow Mach number, Reynolds number, and dynamic pressure to be varied 
independently. Flow conditions used for comparison herein are a Mach number of 0.90, a Reynolds number 
based on the mean aerodynamic chord of 30 x 10 6 , and a dynamic pressure of 1005 psf. The variable 
dynamic pressure enables an extrapolation of the data to a constant dynamic pressure, which eliminates the 
effects of static aeroelastic model deformation from the wind tunnel data. In addition, the dynamic pressure 
step is required to connect low Reynolds number data to high Reynolds number data (10 x 10 6 to 80 x 10 6 ) 
due to facility limitations which prohibit the full Reynolds number sweep to be performed at a constant 
dynamic pressure. 

Data presented were acquired in two separate tests: the first test acquired surface pressure data, and the 
second test acquired force-and-moment data without the potentially adverse influence of pressure 
instrumentation on the balance. Uncertainty of the pressure coefficients was estimated to be less than ±0.01 
at a 95-percent confidence level. The observed variability in the lift, drag, and pitching-moment coefficients 
was ±0.003, ±0.0003, and ±0.0005, respectively. 


Computational Method 


Computational Algorithm 

The computational method used in this study is version 4 of the CFL3D Navier-Stokes code described in 
detail in reference 2. This code, developed at Langley Research Center, solves the three-dimensional 
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time-dependent thin-layer Navier-Stokes equations with a finite-volume formulation. The equations are writ- 
ten in conservation form with the generalized curvilinear coordinates q, and £ as 


where 


dt d$ dT] 


+ ^(H-H V ) = 0 


Q = 



p u pv pw e ] 


( 1 ) 


( 2 ) 


The coordinates q, and £ are oriented essentially in the streamwise, spanwise, and normal directions, 
respectively. It is assumed that the C, coordinate lines are nearly orthogonal to the body at the body surface. 
The symbols t, p, w, v, w, and e denote the time, density, Cartesian velocity components, and total energy per 
unit volume, respectively; F, G, and H are the inviscid fluxes in the £,, q, and £ directions, respectively; and 
H v is the viscous flux in the £ direction. The viscous fluxes, F v and G v , are dropped from the governing equa- 
tions because of the thin-layer assumption, in which the dominant effects are assumed to arise from the vis- 
cous diffusion normal to the body surface. The state equations are written assuming an ideal gas. Molecular 
viscosity is calculated by using Sutherland’s law and Stoke’s hypothesis. The Reynolds stresses are modeled 
using the four turbulence models described in the Turbulence Models section. The governing equations are 
discretized to be consistent with conservation laws in integral form. 

A second-order accurate upwind-biased spatial differencing scheme is used for solving these equations. 
Roe’s flux-difference splitting (ref. 8) is used to construct the upwind differences for the convective and the 
pressure terms. The spatial derivatives are written conservatively as a flux balance across the cell, and the 
shear stress and heat transfer terms are centrally differenced. Spatial approximate factorization (AF) and 
Euler backward integration after linearization in time results in the solution through 5 by 5 block-tridiagonal 
matrix inversions in three directions. An approximate diagonal form of the spatial factors is employed to 
reduce computational time. Convergence acceleration is obtained by using a multigrid full-approximation 
scheme (FAS) and mesh sequencing for flow initialization. 


Turbulence Models 


Baldwin-Lomax with Degani-Schiff modification. The Baldwin-Lomax turbulence model (ref. 4) is 
used widely throughout the CFD community; its capabilities and limitations are well-known. In short, it is 
generally considered a good model for the prediction of attached flows, but it is deficient for flows with any 
significant separated regions. In particular, the Baldwin-Lomax model tends to predict shocks too far down- 
stream for separated transonic flows over aerodynamic configurations. Degani and Schiff (ref. 5) modified 
the original Baldwin-Lomax model to extend the model in a rational manner to permit an accurate determina- 
tion of the viscous length scale for high-angle-of-attack flows in regions of cross-flow separation in which a 
strong leeward vortical flow structure exists. 

The original Baldwin-Lomax model is a two-layer algebraic model. The inner eddy viscosity is deter- 
mined via 


where 


H, ;i = p(0.4Z>>’) 2 £2 


(3) 


( 4 ) 
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and 


+ ^ a/Pw^W/d t \ 

y (Baldwin-Lomax) = [7 ^ 

^ W (B-L) 

where y is the distance to the wall or wake cut, Q is the magnitude of the vorticity, and the subscript w refers 
to a wall value. 

The outer eddy viscosity is given by 

p r>0 = (0.0168)(1.6)pF wake r ( 6 ) 

where 


F wake = 111111 


v F , 

-max max’ 


v /W:«, - lu { U s 

max \^J i fj i i n 


(7) 


and 



( 8 ) 


and 


F = y&D (9) 

In wakes, D is assumed to be 1. The subscript max in equations (7) and (8) refers to the value at the loca- 
tion where F is at its maximum along grid lines that are oriented normal to the body (i.e., in the viscous direc- 
tion for thin-layer Navier-Stokes). 

The eddy viscosity p, is determined by marching away from the wall along these same grid lines and is 
taken as \i t i from the wall (where the value is zero) to the point above the wall where p, , first exceeds p, 0 . 
Thereafter p, is given the value p t 0 . Transition is modeled by setting p, to zero along all grid lines that are 
normal to the wall within a preselected range of laminar flow. 

The major difficulty encountered in applying the Baldwin-Lomax turbulence model to bodies with 
cross-flow separation is that of properly evaluating the scale length y max and, in turn, of determining p, 0 for 
boundary layer profiles in the cross-flow separation region. This difficulty arises because multiple points of 
local maximum can be obtained in vortical flow regions. To eliminate these difficulties, Degani and Schiff 
modified the Baldwin-Lomax turbulence model so that at each axial station the code searches radially out- 
ward along successive rays, sweeping from the windward to the leeward plane of symmetry. Along each ray, 
the code sweeps outward to find the first peak in F(y) and cuts off the search when the peak is reached. To 
prevent the selection of extraneous peaks that might be caused by an unsteady behavior in F(y), a peak is 
considered to have been found when the value of F(y) drops to 90 percent of the local maximum value. 

As proposed by Baldwin and Lomax, the turbulence model examines a quantity containing the local 
fluid vorticity magnitude to determine the length scale and thus the eddy viscosity coefficient. The modifica- 
tions madp by Degani and Schiff (ref. 5) permit the model to differentiate between the vorticity within the 
attached boundary layers from the vorticity on the surfaces of separation and thus to select a length scale 
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based on the thickness of the attached boundary layers rather than one based on the radial distance between 
the body surface and the surface of separation. The latter scale length can be much greater than that based on 
the boundary layer thickness, and if used, the resulting eddy-viscosity coefficient can be as much as two 
orders of magnitude too high. In general, the eddy-viscosity coefficient will cause the primary vortices to be 
smaller than those observed experimentally, and the primary separation line will be located closer to the 
leeward symmetry plane. In addition, the secondary (and tertiary) separation lines will not occur in the com- 
puted flows. 


Baldwin-Barth. The Baldwin-Barth (B-B) turbulence model (ref. 6) is a one-equation turbulence model 
derived from a simplified form of the k-e equations. The model solves a partial differential equation (PDE) 
over the whole field for the turbulent Reynolds number R T , which is directly related to the kinematic eddy 
viscosity v f : 


D\ vRe 7 


! (B-B). 


Dt 
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v + v 

V 0 £ (b-b) J . 
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Then 


|i r = 0.09pfvRe r | D, D 2 . (11) 

‘ V. 1 (B-B) J 1 (B-B) Z (B-B) 

The variables v f , a , f 7 , P m R D, , and D 9 are described in detail in reference 6. In 

t (B-B) Z (B-B) (V-V) ‘(B-B) Z (B-B) 

this formulation of the B-B turbulence model, the thin-layer assumption has been used for the source term 
(the last term in the PDE). 

This equation is solved implicitly by using 3-factor AF, with first-order upwind differencing used on the 
advective terms. Local time stepping is employed to accelerate convergence. The transition location is mod- 
eled by phasing out the source term along all grid lines that are normal to the wall within a preselected range 
of laminar flow. 


Spalart-Allmaras. The Spalart-Allmaras (S-A) turbulence model is a one-equation turbulence model 
derived by “using empiricism and arguments of dimensional analysis, Galilean invariance, and selective 
dependence on the molecular viscosity” (ref. 3). The model solves a PDE over the whole field for a working 
variable V( S A ^ related to the eddy viscosity: 


P ^(S-A) 

Dt 


-0.933v (S . A) V 2 v (S . A)+ 1.5V - [(v + 1 . 622 v (S . a) )Vv (S . a) ] 
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Then 


V-t 


P^(S-A)(W V ) 3 

(v ( s-a / v ) 3 + 357 - 911 


(13) 


The terms , and /„, are described in detail in reference 3. The term rf/c A , denotes 

the nearest distance to any wall. The S-A model is very similar in form to the B-B model, although the S-A 
model includes a destruction term that is not present in the other model. This lack of a destruction term in the 
B-B model is responsible for a mild inconsistency in modeling isotropic turbulence and could also invalidate 
the model in the class of shear flows in which v r decreases (such as an axisymmetric wake). The PDE is 
solved with the same implicit method used in the B-B model, and transition is modeled in a similar fashion as 
well. 


Menter’s SST. The Menter’s Shear Stress Transport (SST) turbulence model is a modified version of the 
k - oo two-equation turbulence model. This model is based on Bradshaw’s assumption that the principal shear- 
stress is proportional to the turbulent kinetic energy, which is introduced into the definition of the eddy vis- 
cosity. TTiis model transforms a k-e model into a k-co formulation, and then the original equation is multi- 
plied by a function (sst> and die transformed model is multiplied by a function (1-F). Both are added 
together, giving the following equations: 
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(14) 


(15) 


where the constants P(SST)’ 7 (SST)’ °^(SST)> ctco (SST)’ P(SST)’ 411(1 ^i (SST) ^ described in detail in 
reference 7. The eddy viscosity is defined as 


Pr = 


pa, k 

1 (SST) 


(16) 


maxf a-, co ;QF~ j 

V a (sst) z (sst)J 

where Q is the magnitude of the vorticity; F 2 (sst> and ai(ssT) are also described in detail in reference 7. 


Grid Generation 

The grid generated about the NTF 2.2-percent HSR Reference H baseline geometry (fig. 1) was created 
using the GRIDGEN grid generation package (ref. 9). Three grids were generated for this geometry in order 
to assess the effect of grid topology and the representation of the actual wind tunnel model geometry. The 
different topologies and geometric representations were chosen because these topologies and geometric rep- 
resentations are the most common ones used for this configuration. 
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The first grid was a two-block grid with an O-C (O in the streamwise direction and C in the spanwise 
direction) topology with a 5-in. full-scale transition from blunt-to-sharp leading edge at the crank location on 
the wing (fig. 4). This grid had 105 points in the streamwise direction, 177 points in the spanwise direction, 
and 81 points normal to the surface in the first block, which defines the wing-body configuration; the second 
block, which defines the sting, had 49 points in the streamwise direction, 177 points in the spanwise 
direction, and 81 points normal to the surface. The total number of grid points was 2207898. The normal 
spacing adjacent to the surface was 5 x 10 -5 over the entire surface. The surface spacing distribution corre- 
sponds to a nondimensional y + value of approximately 1.6 for M = 0.90 and Re a = 30 x 10 6 . The far-field 
boundary extends to 20c in the circumferential and upstream directions, and the second block extends 
approximately 20c downstream of the first block. 

The second grid used was a single block C-0 (C in the streamwise direction and O in the spanwise direc- 
tion) topology grid with a 5-in. full-scale transition from blunt-to-sharp wing leading edge at the crank; this 
grid defines both the wing-body configuration and the sting. This grid had 133 points in the spanwise direc- 
tion, 249 points in the streamwise direction, and 81 points in the normal direction. The total number of grid 
points was 2682477. The normal spacing was 5 x 10 -5 over the entire surface. The surface spacing distribu- 
tion corresponds to a nondimensional y + value of approximately 1.0 for M= 0.90 and Re- = 30 x 10 -6 . The 
far-field boundary extends to 20c in the circumferential, upstream, and downstream (of the wing-body con- 
figuration) directions. Figure 4 shows the differences in the surface grids of the O-C grid on the wing with a 
5-in. full-scale transition from blunt-to-sharp wing leading edge and the C-0 grid on the wing with a 5-in. 
full-scale transition from blunt-to-sharp wing leading edge. Figure 5 shows a representation of the volume 
grid for the C-0 grid topology with a 1-in. full-scale transition from blunt-to-sharp wing leading edge. Fig- 
ures 6 and 7 show that both grids are converged for lift and that the C-0 grid is converged for drag at the 
given conditions. The poor convergence quality of the O-C grid is believed to be caused by a lack in the grid 
quality. The results presented in this paper are for the finest mesh, and the differences due to the mesh topol- 
ogy will be discussed in the Results and Discussion section. 

The wind tunnel model had a discontinuity between the sharp wing leading edge and the blunt wing 
leading edge at the crank. In an attempt to model this discontinuity, a third grid was generated that was iden- 
tical to the second grid (C-0 topology) with the exception of a 1-in., rather than a 5-in., full-scale geometric 
transition region between the blunt and sharp wing leading edges at the crank. Figure 8 shows the differences 
in the surface grids of the C-0 grid on the wing with a 5-in. full-scale transition from blunt-to-sharp wing 
leading edge and the C-0 grid on the wing with a 1-in. full-scale transition from blunt-to-sharp wing leading 
edge. The sketch shows where the geometric wedge is on the wing and what is meant by a 1-in. versus a 5-in. 
full-scale geometric transition from blunt-to-sharp wing leading edge. 

Results and Discussion 

Grid Topology and Wedge Size 

A grid study was performed with the three grids described previously to determine which grid gave the 
best comparisons with the experimental data. All runs made for this grid study were performed at 
M = 0.90, Re- = 30 x 10 6 , arid at a = 5°. These runs were all made with CFL3D by using the S-A turbu- 
lence model. Figure 9 shows the pressure distributions at the six spanwise stations on the wing of the HSR 
2.2-percent NTF geometry for the C-0 topology, the O-C topology, and the experimental data. Figure 10 
shows the pressure contours for the two topologies. As shown in figure 9, both topologies show similar 
agreement with the experimental data, except aft of the crank of the wing (x = 46.5 in.) where the C-0 topol- 
ogy picks up more of the characteristics of the experimental pressure distribution than does the O-C topol- 
ogy. The pressure contours in figure 10 show that the C-0 grid topology yields a different surface pressure 
distribution than the O-C grid topology, primarily on the outboard wing panel. 
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Figure 1 1 shows the pressure distributions at the six spanwise stations on the wing of the HSR 
2.2-percent NTF geometry for the 1-in. full-scale wedge versus the 5-in. full-scale wedge. Figure 12 shows 
the pressure contours of the two wedge sizes. As shown in figures 11 and 12, the wedge size does not have a 
significant effect on the pressure distributions or the pressure contours. After we looked at the results from 
these two studies, the single block C-0 grid with the 1-in. wedge was used for the turbulence model study 
because the 1-in. wedge more closely represents the actual wind tunnel model geometry. 

Code Modifications 

During this study, the CFL3D Navier-Stokes code was modified to increase the calculation accuracy of 
the minimum distance function by incorporating a surface interpolation scheme based on normal projection. 
The prior approach used the nearest body cell center point. This modification only affected the results of the 
S-A and SST turbulence models. The minimum distance function affects the L, and d ( s _a^ terms in 

W ( S - A ) '■‘ 3 

equation (12) for the S-A model and the F 1(SST) term in equation (15) for the SST model. Figures 13 and 14 
show that the modification made to the code does not significantly affect the pressure distributions when 
either the S-A turbulence model or Menter’s SST turbulence model is used, respectively. Figure 15, how- 
ever, indicates that the forces are affected by this modification to CFL3D. This change is seen because the 
modification changes the solution, which, in turn, results in a change in the integrated forces. For the S-A 
case, the lift is increased by approximately 0.004 and the drag is decreased by approximately 20 drag counts. 
For Menter’s SST results, the lift remains about the same, but the drag is increased by 3 drag counts. The 
change to the code brought the two calculations into closer agreement with each other but took the S-A 
results farther away from the experimental data. 

This modification also decreased the amount of computer time needed to obtain a converged solution by 
approximately 40 percent for the S-A case. The SST case was restarted from the previous runs with the new 
code, while the S-A case was completely reran with the new code; therefore, the amount of time savings for 
the SST case is unknown, but it is expected to decrease the amount of computational time by about the same 
amount as observed with the S-A model. 

Turbulence Model Comparison 

Four different turbulence models, the equilibrium model of Baldwin-Lomax with the Degani-Schiff 
modifications, the one-equation Baldwin-Barth model, the one-equation Spalart-Allmaras model, and 
Menter’s two-equation Shear Stress Transport model were used to perform runs for the HSR 2.2-percent 
NTF geometry at M = 0.90, Re s = 30 x 10 6 , and at three angles of attack, 1°, 5°, and 10°. In this investiga- 
tion, it is assumed that the entire boundary layer flow is turbulent. All computations were performed with 
CFL3D using a C-0 topology grid with approximately 2.7 million grid points. 

The solutions were considered to be converged when the drag coefficient changed less than one drag 
count over 100 iterations and when the residual was decreased by at least three orders of magnitude. Typical 
convergence histories for all three angles of attack are shown in figures 16 through 18; the run procedure, 
along with the ran time, is given for the residual history for each case. Fust, the flow is initialized on two 
coarser grids, and then the solution is interpolated onto the finest level mesh where iterations are performed 
until convergence is reached. Multigrid iterations are implemented on the level 2 (medium coarse) mesh and 
the level 3 (fine) mesh in order to increase the convergence speed. The a = 1.0° case shown is with the SST 
model, the a = 5.0° case shown is with the B-B model, and the a = 10.0° case shown is with the S-A model. 
The Baldwin-Lomax with Degani-Schiff (B-L) runs were made on the Cray Y-MP at Langley Research 
Center, while the runs for the other three turbulence models were performed on the Numerical Aerospace 
Simulator (NAS) Cray C-90. All run times given in the figures are in C-90 equivalent central processing unit 
(CPU) hours. The jump seen for the last few iterations in the drag histories in figures 16 and 17 was caused 
by a correction to the method of calculating the drag coefficient in CFL3D during the course of this study. 
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This correction involved a change in the way the forces were integrated; therefore, the flow field itself was 
not affected. Since these cases were already converged when the correction was made, a few (5 to 10) itera- 
tions were run on each case, even though only one iteration was needed in order to obtain the correct drag 
values for these cases. 

Table 1 gives a tabulated listing of the run times for all four turbulence models at all three angles of 
attack; all run times given are in C-90 equivalent CPU hours. The results in the table show that S-A uses the 
least ps/gridpoint/iteration at all three angles of attack, with B-L using the next fewest at all three angles of 
attack. The results in the table also show that S-A converges the fastest, in terms of CPU hour usage, of the 
four models at all three angles of attack. At a = 1.0° and a = 10.0°, the SST model has the next fastest con- 
vergence rate, but at a = 5.0° B-L has the second fastest convergence rate. Finally, the results in the table 
show that S-A uses the lowest percentage of CPU hours on the fine mesh of the four models at all three 
angles of attack. At a = 1.0° and a = 10.0°, the B-B model uses the next to the lowest percentage of CPU 
hours on the fine mesh, and at a = 5.0°, B-L uses the second lowest percentage of CPU hours on the fine 
mesh. These percentages are important when considering that B-L needs 120 MW of memory to run on the 
fine mesh, while S-A, SST, and B-B require 127 MW, 127 MW, and 136 MW, respectively. The amount of 
memory each model uses is important because of the limited availability of computer resources. 

Pressure comparisons. The spanwise and chordwise pressure distributions for the a = 1.0° case are 
shown in figures 19 and 20, respectively. These figures indicate that for this low angle of attack, all turbu- 
lence models agreed well with the experimental data at all six spanwise and at all three chordwise pressure 
stations. This behavior is expected because this relatively benign flow is mostly attached, and all models 
should behave the same for attached flow. The spanwise and chordwise pressure distributions for the 
a = 5.0° case are shown in figures 21 and 22, respectively. The spanwise pressure distributions for the 
cl = 5.0° case shown in figure 21 indicate that up to the crank in the wing, all turbulence models predict the 
flow characteristics very well. At the crank region (x = 46.5 in.) none of the models pick up the detail seen in 
the outboard 20-percent span in the experimental data. At an x station of 49.55 in., B-L is the only model that 
picks up the apparent suction peak at 90-percent span in the experimental data. At the x = 53.16-in. station, 
all models perform equally well, with B-L showing a suction peak at approximately 80-percent span. The 
experimental pressures do not extend this far outboard; therefore, it is difficult to determine whether the suc- 
tion peak is real. Figure 22 shows that for the a = 5.0° case, all four models give good agreement with exper- 
imental data at the first two y locations. At the y = 9.95-in. location, no models match the pressure 
measurement on the upper surface. 

The spanwise and chordwise pressure distributions for the a = 10.0° case are shown in figures 23 and 24, 
respectively. The spanwise pressure distributions for the a = 10.0° case shown in figure 23 indicate that at 
the x station of 34.245, all models predict the flow similarly, and this prediction varies slightly from the 
experimental flow. At the next two stations, that is, at x = 37.98 in. and x = 41.45 in., no models pick up the 
magnitude of the suction peak seen at approximately 90-percent span and 85-percent span, respectively, in 
the experimental flow. At the crank region, that is, at x = 46.5 in., the B-L model overpredicts the suction 
peak seen at approximately 80-percent span in the experimental data by approximately the same amount that 
the other three models underpredict the height of the suction peak. At an x station of 49.55 in., B-L comes the 
closest to matching the fluctuations in the experimental data. At the x = 53.16-in. station, B-L comes the clos- 
est of the four models to picking up the upper surface pressure distribution seen in the experiment. Figure 24 
shows that for the a = 10.0° case, no models pick up the dip in the pressure curve which occurs at 60-percent 
chord in the experimental data at the y = 3.036-in. station. At the second y location, y = 5.544 in., B-L comes 
the closest to picking up the drop in the pressure curve at 30-percent chord, but all four models pick up the 
pressure at 60-percent chord. At the y = 9.95-in. location, all models come close to predicting the experimen- 
tal value for the upper surface pressure; B-L shows some fluctuations in the pressure distribution which the 
other models do not predict. There is no way of knowing whether these fluctuations really do occur in nature 
because there were no pressure values taken experimentally. 
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Surface streamlines. Figures 25 through 27 show the computational surface streamlines for the a = 1 .0°, 
5.0°, and 10.0° cases, respectively, for all four turbulence models. In these figures, a separation line is indi- 
cated by the streamlines coalescing to a line, while an attachment line is a line from which the streamlines are 
flowing. As shown in figure 25, the computational surface streamlines are aligned with the free-stream flow 
on both the outboard and inboard wing sections for all four turbulence models at a = 1.0°. The data indicate 
there is no separation present at this low angle of attack. In figure 26, all computational surface streamlines 
for the a = 5.0° case show a separation at the crank with an attachment near the outboard wing leading edge. 
The B-L solution shows an additional separation and reattachment at the outboard wing leading edge that the 
other models do not show. Figure 27 shows the computational streamlines for the a = 10.0° case. All models 
show a separation at the crank region of the wing. The B-L case, however, also shows a separation at the 
inboard leading edge of the wing and a reattachment line in the outboard region of the wing. The other three 
models show a separation region at the tip of the wing that the B-L model does not show. No experimental 
data are available for comparison. 

Force and moment comparisons. Figures 28 through 31 show the drag polars, lift-to-drag ratios, lift 
curves, and pitching-moment curves, respectively, for the four different turbulence models and the experi- 
mental data. Table 2 shows the data which are plotted, and table 3 gives a numerical breakdown of the pres- 
sure drag and viscous drag components for all calculations made for this study. Figure 28 shows the drag 
polar for the four different turbulence models and the experimental data. At ct = 1.0°, all turbulence models 
predict the lift-to-drag ratio to within the accuracy of the experimental drag value (i.e., ±0.0003). At the 
design point, a = 5.0°, B-B is approximately 8 drag counts low, with S-A, SST, and B-L being 10, 11, and 
15 drag counts low, respectively. At a = 10.0°, B-L predicts lift and drag lower than the other three models, 
but it does come the closest to predicting the ratio of lift to drag. 

The c omp utational wing pressure distributions discussed earlier in this paper indicated that B-L agreed 
the best with the experimental data while the drag data in figure 28 indicate that the B-B model agrees some- 
what better with the experiment at a = 5.0°. This discrepancy is probably caused by a difference in the vis- 
cous drag components of the turbulence model results, which are shown in table 3. The B-L turbulence 
model viscous-drag-component results are approximately 10 drag counts lower than the other three models at 
all three angles of attack. At a = 10.0°, the pressure component of drag for B-L is approximately 60 drag 
counts lower than the other three models. This factor can account for the results seen at this angle of attack. 
Note also that the equivalent flat plate viscous drag = 0.0064 for M — 0.90 and Re- = 30 x 10'\ which is 
7 to 17 counts higher than the computations herein. Figure 32 shows that if the flat plate equivalent 
viscous-drag component is substituted for the computational viscous-drag component, the computational 
data agree much better with the experimental data. 

The lift-to-drag ratios for the four turbulence models versus the experimental data are shown in fig- 
ure 29. This figure shows that at ct = 1.0°, all four models predict approximately the same L/D value, which 
agrees well with the experimental value. At a = 5.0°, B-B comes closest to matching the experimental data. 
At a = 10.0°, B-L predicts the experimental L/D level closer than the other three models. The discrepancies 
between the computational and experimental pressure data previously examined and the force data compari- 
sons shown here are consistent with the differences in the viscous drag and pressure drag components of the 
four turbulence models. 

Figure 30 shows the lift curve for the four turbulence models and the experimental data. At a = 1.0°, all 
models are slightly higher than the experimental data. At the design point, ct — 5.0°, B-L is closer to the 
experimental lift curve than the other three models, and at ct = 10.0°, B-L is again closer to the experimental 
lift curve than the other three models. 

Figure 31 shows the pitching-moment curves for the four turbulence models and the experimental data. 
At a = 1.0°, all models underpredict the pitching moment. At a = 5.0°, all models again underpredict the 


12 



pitching moment, but B-L comes the closest to the experimental values. At a = 10.0°, all models again 
underpredict the pitching moment, but B-L comes the closest to the experimental values. 

These results show that the computed forces and pressure compared reasonably well with the experi- 
mental data, with the B-L and B-B models showing superior agreement. None of the turbulence models pre- 
dict the experimental flow well at all conditions tested herein. B-L seems to give the best agreement with the 
force data at the higher angles of attack, while B-B agrees better with the force data at lower angles of attack. 
Also, substitution of the theoretical equivalent flat plate C Dv for the Navier-Stokes computed signifi- 
cantly improves the agreement between the computational and the experimental data. 

Conclusions 

Four turbulence models, the equilibrium model of Baldwin-Lomax (B-L) with the Degani-Schiff (D-S) 
modifications, the one-equation Baldwin-Barth (B-B) model, the one-equation Spalart-Allmaras (S-A) 
model, and Menter’s two-equation Shear Stress Transport model (SST) have been evaluated for transonic 
flow over a high-speed research (HSR) geometry using the thin-layer, upwind, Navier-Stokes flow solver 
known as CFL3D. The flow conditions tested were a Mach number of 0.90 and a Reynolds number of 
30 x 10 6 based on chord for a range of angles of attack of 1°, 5°, and 10° which correspond to tests per- 
formed in the National Transonic Facility (NTF) at Langley Research Center. The results indicate that either 
a C-O topology or an O-C topology will give reasonable results for this configuration, but the C-0 topology 
was found to be slightly better. In addition, the detail of the geometric modeling of the blunt-to-sharp 
leading-edge transition at the wing crank did not significantly affect the results, at least for the variation of 
detail examined herein. A distance function modification within CFL3D improved correlation between the 
S-A and Menter’s SST turbulence model results. The turbulence model study showed that the computed 
forces and pressure compared reasonably well with the experimental data. Finally, the results also showed 
that S-A converges the fastest of the four models at all three angles of attack. 


NASA Langley Research Center 

Hampton, VA 23681-2199 

June 17, 1997 
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Table 1. Tabulated Run Times 







CPU hr used on 



a. 

Turbulence 

ps/gridpoint/ 

CPU hr 

fine mesh, 



deg 

model 

iteration 

(C-90 equiv.) 

percent 



1 

B-L 

13.12 

6.77 

87 




1 

B-B 

13.99 

7.80 

80 




1 

S-A 

10.14 

4.90 

77 




1 

SST 

14.54 

6.29 

86 




5 

B-L 

11.05 

5.52 

89 




5 

B-B 

12.12 

9.74 

93 




5 

S-A 

10.18 

4.91 

77 




5 

SST 

12.56 

11.31 

91 




10 

B-L 

10.75 

10.58 

91 




10 

B-B 

14.10 

10.64 

89 




10 

S-A 

10.19 

4.92 

77 




10 

SST 

13.97 

10.29 

91 




Table 2. Tabulated Values of Drag, Lift, and Pitching Moment 




a, 

Exp. 

Baldwin- 

Baldwin- 

Spalart- 

Menter’s 

Force 

deg 

II 

>— k 

8 

Lf \ 

Lomax 

Barth 

Allmaras 

SST 


1 

0.0074 

0.0054 

0.0062 

0.0063 

0.0060 

C D 


5 

.0183 

.0176 

.0191 

.0188 


.0185 


10 

.0751 

.0768 

.0847 

.0855 


.0855 


1 

-.008 

.001 

.002 

.002 


.001 

C, 

5 

.217 

.226 

.232 

.232 


.230 


10 

.494 

.496 

.542 

.548 


.550 


1 

.014 

.011 

.011 

.011 


.011 

r ... 

5 

.008 

.006 

.004 

.004 


.004 


10 

.006 

.004 

-.007 

-.009 


.010 




Table 3. Tabulated Values of Computational Pressure and Viscous Drag 
[Equivalent Flat-Plate Viscous Drag; M = 0.90; Re a = 30 x 10 6 ; C Dv = 0.0064] 


Force 

a, 

deg 

Baldwin- 

Lomax 

Baldwin- 

Barth 

Spalart- 

Allmaras 

Menter’s 

SST 


1 

0.0004 

0.0005 

0.0005 

0.0005 

c Dp 

5 

.0130 

.0135 

.0133 

.0131 


10 

.0723 

.0789 

.0800 

.0801 


1 

.0049 

.0057 

.0058 

.0054 

c Dv 

5 

.0046 

.0056 

.0055 

.0053 


10 

.0045 

.0057 

.0055 

.0052 
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Figure 4. Wire-frame surface grids of O-C and C-0 topologies. 
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Figure 5. Representation of the C-0 type volume grid for NTF 2.2-percent HSR geometry. 


20 







1 /Total number of points 


Figure 7. Grid density study showing computed lift coefficient plotted against inverse of mesh size for 
O-C and C-0 topologies. M = 0.90; a = 5.0°; Reg = 30 x 10' . 
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ylb/2 


y/b/2 


• Exp., q = 1005 psf, upper surface 

O Exp., q = 1005 psf, lower surface 

O-C topology, upper surface 

— — — O-C topology, lower surface 

C-0 topology, upper surface 

C-0 topology, lower surface 

Figure 9. Surface pressure comparison between experiment, O-C, and C-0 grid topology results. 
M = 0.90; a = 5.0°; Reg = 30 x 10 6 . Dimensions are in inches. 
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O-C topology 



Figure 10. Pressure contour comparison between O-C and C-0 grid topologies for Spalart-Allmaras 
turbulence model. M= 0.90; a = 5.0°; Re^ = 30 x 10 6 . 
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ylbl2 ylbl 2 

• Exp., Q = 1005 psf, upper surface 

O Exp., q = 1005 psf, lower surface 

1-in. wedge, upper surface 

— — — 1 -in. wedge, lower surface 

5-in. wedge, upper surface 

5-in. wedge, lower surface 

Figure 11. Surface pressure comparison between experiment, 1-in., and 5-in. full-scale wedge results. 
M = 0.90; a = 5.0°; Reg = 30 x 10 6 . Dimensions are in inches. 
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turbulence model. M = 0.90; a = 5.0°; Reg = 30 x 10 6 . 
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• Exp., q = 1005 psf, upper surface 

O Exp., q = 1005 psf, lower surface 

Original CFL3D, S-A model, upper surface 

— — — Original CFL3D, S-A model, lower surface 

Modified CFL3D, S-A model, upper surface 

Modified CFL3D, S-A model, lower surface 

Figure 13. Code comparison, with and without modified distance function, for Spalart-Allmaras 
turbulence model. M = 0.90; a = 5.0°; Re a = 30 x 10 6 . Dimensions are in inches. 
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y /b/2 ylb/2 

• Exp., q = 1005 psf, upper surface 

O Exp., q = 1005 psf, lower surface 

Original CFL3D, SST model, upper surface 

— — — Original CFL3D, SST model, lower surface 

Modified CFL3D SST model, upper surface 

Modified CFL3D SST model, lower surface 

Figure 14. Code comparison, with and without modified distance function for Menter’s SST turbulence 
model. M = 0.90; a = 5.0°; Re? = 30 x 10 6 . Dimensions are in inches. 
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(a) Residua] history. 

Figure 16. Typical convergence history for a =1.0° case, SST turbulence model; total CPU 
hours = 6.29. 


31 


.0035 


.0034 



000 l i_i 1 1 1 i 1 1 1 1 i 1 1 1 i 1 i 1 1 1 

0 500 1000 1500 2000 

Iterations 


(b) Drag history. 
Figure 16. Concluded. 
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(a) Residual history. 

Figure 17. Typical convergence history for a = 5.0° case, B-B turbulence model; total CPU 
hours = 9.74. 
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0.04 CPU hr 
single-cycle 
level 1 mesh 


1.09 CPU hr 
multi-cycle 
level 2 mesh 


3.8 CPU hr 
multi-cycle 
level 3 mesh 



(a) Residual history. 

Figure 18. Typical convergence history for a = 10.0° case, S-A turbulence model; total CPU 
hours = 4.9. 
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(b) Drag history. 
Figure 18. Concluded 
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Figure 19. Spanwise wing pressure distributions. M = 0.90; a = 1.0°; Re s = 30 x 10 6 . Dimensions are 
in inches. 
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Figure 20. Chordwise wing pressure distributions. M- 0.90; a = 1.0' 
in inches. 
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; Re e = 30 x 10 6 . Dimensions are 
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ylb/2 ylbll 


• Exp., q = 1005 psf, upper surface 

O Exp., q = 1005 psf, lower surface 

Spalart-Allmaras, upper surface 

— — — Spalart-Allmaras, lower surface 

Baldwin-Barth, upper surface 

Baldwin-Barth, lower surface 
Baldwin-Lomax, upper surface 

— — — — Baldwin-Lomax, lower surface 

SST, upper surface 
- - SST, lower surface 

Figure 21. Spanwise wing pressure distributions. M = 0.90; a = 5.0°; Re s = 30 x 10 6 . Dimensions are 
in inches. 
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Figure 22. Chordwise wing pressure distributions. M = 0.90; a = 5.0° 
in inches. 
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; Re^ = 30 x 10 6 . Dimensions are 
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y/b/2 




• Exp., q = 1005 psf, upper surface 

O Exp., q = 1005 psf, lower surface 

Spalart- Allmaras, upper surface 

— — — Spalart- Allmaras, lower surface 

Baldwin -Barth, upper surface 

Baldwin-Barth, lower surface 
Baldwin-Lomax, upper surface 

— — — — Baldwin-Lomax, lower surface 

SST, upper surface 
SST, lower surface 

Figure 23. Spanwise wing pressure distributions. M = 0.90; a = 10.0°; Re^ = 30 x 10 6 . Dimensions are 
in inches. 
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Figure 24. Chordwise wing pressure distributions, 
in inches. 


42 


• Exp., q = 1005 psf, upper surface 

O Exp., q = 1005 psf, lower surface 

Spalart-Allmaras, upper surface 

— — — Spalart-Allmaras, lower surface 

Baldwin-Barth, upper surface 

Baldwin-Barth, lower surface 
Baldwin-Lomax, upper surface 

— — — — Baldwin-Lomax, lower surface 


SST, upper surface 
- - SST, lower surface 


= 0.90; a = 5.0°; Re 5 = 30 x 10 6 . Dimensions are 
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Figure 31. Pitching- moment 
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Figure 32. Drag polar using flat-plate equivalent viscous drag instead of computational viscous drag 
component. M = 0.90; Reg = 30 x 10 6 . 
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